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Orbital evolution of an interplanetary dust particle under action of an interstellar gas flow is investigated. Secular time 
derivatives of the particle orbital elements, for arbitrary orbit orientation, are presented. An important result concerns 
secular evolution of semi- major axis. Secular semi- major axis of the particle on a bound orbit decreases under the action 
of fast interstellar gas flow. Possible types of evolution of other Keplerian orbital elements are discussed. The paper 
compares influences of the Poynting-Robertson effect, the radial solar wind and the interstellar gas flow on dynamics 
of the dust particle in outer planetary region of the Solar System and beyond it, up to 100 AU. 

Evolution of putative dust ring in the zone of the Edgeworth-Kuiper belt is studied. Also non-radial solar wind and 
gravitational effect of major planets may play an important role. Low inclination orbits of micron-sized dust particles in 
the belt are not stable due to fast increase of eccentricity caused by the interstellar gas flow and subsequent planetary 
perturbations - the increase of eccentricity leads to planet crossing orbits of the particles. 

Gravitational and non-gravitational effects are treated in a way which fully respects physics. As a consequence, some 
of the published results turned out to be incorrect. Moreover, the paper treats the problem in a more general way than 
it has been presented up to now. 

The influence of the fast interstellar neutral gas flow might not be ignored in modeling of evolution of dust particles 
beyond planets. 

Key words. ISM: general, Celestial mechanics, Interplanetary medium 



1. Introduction 

Motion of stars relative to their local interstellar medium 
is frequent /usual process in galaxies. Neutral atoms pene- 
trate into the Solar System due to the relative motion of the 
Sun with respect to the interstellar medium. This flow of 
neutral atoms through a heliosphere has been investigated 
in many papers, e.g. Fahr (1996), Lee et al. (2009), Mobius 
et al. (2009). Motion of dust in interplanetary space can be 
affected by the neutral gas penetrating into the heliosphere. 
Influence of this effect on dynamics of dust particles is usu- 
ally ignored in literature. The Poynting-Robertson effect, 
the radial solar wind and the gravitational perturbation 
of planet(s) are usually taken into account (Sidlichovsky 
& Nesvorny 1994; Liou & Zook 1997; Liou & Zook 1999; 
Kuchner & Holman 2003). 

Scherer (2000) has calculated secular time derivatives 
of angular momentum and Laplace-Runge-Lenz vector of a 
dust particle under the action of interstellar gas flow. But 
Scherer's calculations contain several incorrectnesses. He 
has come to the conclusion that semi-major axis of the dust 
particle increases exponentially (Scherer 2000, p. 334). This 
paper presents that semi-major axis of the dust particle 
decreases under the action of interstellar gas flow, in the 
framework of the perturbation theory. 

Motion of dust particles in the zone of the Edgeworth- 
Kuiper belt under the action of the interstellar flow of gas 
has been investigated by Klacka et al. (2009a) . The authors 
have calculated secular time derivatives of orbital elements 
only for the case when interstellar gas velocity vector lies 
in the orbital plane of the dust particle and direction of the 



velocity vector is parallel with y-axis. This paper overcomes 
these restrictions. Moreover, it presents some main proper- 
ties of dust dynamics under the action of the interstellar 
gas. 

2. Secular evolution 

Acceleration of a spherical dust particle caused by the flow 
of neutral gas can be given in the form (Scherer 2000) 



^ = - c D 7jj \v - vh\ {v- v h ) 



(1) 



where vh is velocity of the neutral hydrogen atom, v is 
velocity of the dust grain, crj is the drag coefficient, is 
the collision parameter. For the collision parameter we can 
write 



Jh = n H A , 

m 



(2) 



where nin is mass of the neutral hydrogen atom, nn is 
the concentration of interstellar neutral hydrogen atoms, 
A = ttR 2 is the geometrical cross section of the spherical 
dust grain of radius R and mass m. The concentration of 
interstellar hydrogen rin is not constant in the entire helio- 
sphere. For heliocentric distances r less than 4 AU tih de- 
creases precipitously from its value in the outer heliosphere 
toward the Sun, due to ionization (Lee et al. 2009). But in 
the outer heliosphere, r € (30 AU, 80 AU), we can assume 
that the concentration of the neutral hydrogen atoms is 
constant nn — 0.05 cm -3 (Fahr 1996). The same assump- 
tion can be used also behind the solar wind termination 
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shock. The shock was crossed by Voyager 1 at a heliocen- 
tric distance 94 AU and by Voyager 2 at 84 AU (Richardson 
et al. 2008). 

We will assume that the speed of interstellar gas is much 
greater than the speed of the dust grain in the stationary 
frame associated with the Sun (v <C Vh). This approxi- 
mation leads to approximately constant value of cd ~ 2.6 
(Baines et al. 1965; Banaszkiewicz et al. 1994; Klacka et al. 
2009a). 

We want to find influence of the flow of interstellar gas 
on secular evolution of particle's orbit. We will assume that 
the dust particle is under the action of solar gravitation and 
the flow of neutral gas. Hence we have equation of motion 



dv 



-77 = »r - c D -f H v 

at r 6 



vh\ 



vh) ; 



(3) 



where fi — GMq, G is the gravitational constant, Mq is 
mass of the Sun, r is position vector of the dust particle 
with respect to the Sun and r = \r\. The acceleration caused 
by the flow of interstellar gas will be considered as a per- 
turbation acceleration to the central acceleration caused by 
the solar gravity. In order to compute secular time deriva- 
tives of Keplerian orbital elements (a - semi-major axis, 
e - eccentricity, ui - argument of perihelion, S7 - longitude 
of ascending node, i - inclination) we want to use Gauss 
perturbation equations of celestial mechanics. To do this, 
we need to know radial, transversal and normal compo- 
nents of acceleration given by Eq. (1). Orthogonal radial, 
transversal and normal unit vectors of the dust particle on 
the Keplerian orbit are (see Fig. 1 and e.g. Pastor 2009) 



e R 



ex = 



e/v 



(cos fl cos(/ + ui) — sin SI sin(/ + u) cosz , 
sinfi cos(/ + ui) + cos SI sin(/ + u;) cosi , 
sin(/ + cu) sini) , 

(—cos SI sin(/ + w) — sinf2 cos(/ + cj) cosz 
— sin SI sin(/ + ui) + cos S7 cos(/ + ui) cosz , 
cos(/ + ui) sini) , 
(sinfi sin 2, — cosf2 sini, cosi) , 



(4) 



(5) 
(6) 



where / is true anomaly. Thus, we need to calculate the 
values of a R = dv/dt ■ e R , ax = dv/dt ■ ex and a at = 
dv/dt- e^. Velocity of the particle in an elliptical orbit can 
be calculated from 



dr d 

~db ~ ~dt 
e sin/ 



{res) 
df 



where 



1 + e cos / dt 



P 



e R + r e T 



df_ 

dt 



(7) 



1 + e cos / 



and p = a{l 
Kepler's law df /dt 
easily verify that 



2 ). 



(v - v H ) ■ e R 
{v - v H ) ■ e T 



(v - v H ) ■ e 



N 



(8) 

In this calculation also the second 
^TjMp/r 2 must be used. Now, one can 



(9) 
(10) 

(11) 



vh c e sin / — vh • eR. 
vh & e sin / — A , 
vh er (1 + e cos /) — vh • er 
vh a (1 + e cos /) — B , 



vh ■ e N 



C 




Fig. 1. A particle on an elliptical orbit depicted together 
with the radial, transversal and normal unit vectors. Angles 
characterizing position of the particle on the orbit are also 
shown. 



where 



a 



and vh = 
\v - v H \ 2 



(12) 



\vh\- Hence 

2e cos / 



£ {1 
p 



B + e{A sin / + B cos /)] + v 2 H , (13) 



where the identity \/ A 2 + B 2 + C 2 = vh was used. 

If we denote components of the velocity vector of hydro- 
gen gas in the stationary Cartesian frame associated with 
the Sun as vh = (vhx, vhy, vhz), then we obtain 

A sin f + B cos/ = (— cos £7 sinw 

— sinS7 cos a; cosi) Vhx 

+ (— sinfi sin a; 

+ cos Q COS U! cos i) VHY 

+ cosw sini vhz = I ■ (14) 

Now we consider only such orbits for which 

a<l, (15) 

or, more precisely, the value a 2 is negligible in comparison 
with a. This is reasonable for orbits with not very large 
eccentricities, since v <C Vh- Using this approximation, Eqs. 
(13)-(14) yield 



\v-v H \ =v H [l {B + el)] 

v H 



(16) 



For radial, transversal and normal components of accelera- 
tion we obtain from Eq. (1), Eqs. (9)-(ll) and Eq. (16) 



cd Jh v H 



A_ f uel 
vh V v h 
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e sin / + — 2- 



H 



(17) 



(It 



CD JH V H 



B (ael _ l 



B 2 

1 + e cos / H — x- 

V 2 H 



B 



(18) 
(19) 



a N = - c D 1h v H C 1 

y vh 

Now we can use Gauss perturbation equations of celestial 
mechanics to compute time derivatives of orbital elements. 
The perturbation equations have the form 



IP 



da 


2 


~dt ~ 


1 - 


de 


ft 


~dl ~ 


V m 


du 


1 


~dt ~ 


e 



— J— [a R e sin / + a T (1 + e cos/)] 



a R sin / + a,T cos / 



cos/ 



1 + e cos / 



2 + ecos/ 
a fl cos / - a T — sin / 



1 + e cos / 



m p 



dn 

~dt 

di 
dt 



/HP 
r 

fJTp 



a N 



sin(/ + uj) 

(2 TV ; COS? 

sin? 

sin(/ + uS) 



sin! 

a^r cos(/ + w) . 



(20) 



Time average of any quantity g during one orbital period 
T can be computed using 



1 f X 



Q dt — 7, 77. 



271 



27ra 3 / 2 
1 



27raVl - e 2 



r 

271 



.9 



4f 



gr 2 df 



(21) 



2„3 _ 



where the second (y/JIp = r df /dt) and the third (4n a 
\iT 2 ) Kepler's laws were used. From Eqs. (17)-(21) we fi- 
nally obtain for the secular time derivatives of the Keplerian 
orbital elements 



da\ 
dt/ 



2ac D jHV Hx l-(y\l + -Y- 



/de 
\di 



dco 
~dt 




(22) 



(23) 



+ C 



cos / 



3e sin uj 
1-e 2 



■ (S cos uj — I sin uj) 



dn 

~dt 



VH 

cp 7ff vh C 
2sini 



(24) 



3e sinw 



1 



H (S cos uj - 

VH 

CD 1H V H C 



J sin uj) 



(25) 



3e cosoj 



1 



-(S sin a; + Icosw) 



(26) 



where the quantities 

5 = (cos 17 cos w — sin 17 sinw cosi) vhx 

+ (sin 17 cos w + cos 17 sinw cosi) vhy 

+ sin u> sini vhz > 
I = (— cos 17 sinw — sin 17 cos a; cosi) 

+ (—sin 17 sin w + cos 17 cos a; cosi) vhy 

+ cos a; sini d^z , 
C = sin 17 sin i vhx — cos 17 sin i vhy + cos i vhz , (27) 

are values of A = vh ■ e R , B = vh ■ &t and C — vh ■ &n 
at perihelion of particle's orbit (/ = 0), respectively. The 
value of C is a constant on a given oscular orbit. The values 
of S, I and C are depicted in Fig. 2. One can use also the 




Fig. 2. A schematic representation of the values S, I and 
C for a given orbit. 



relations 

S cos uj — I sin uj 
S sin uj + I cos uj 



epA • epE 
e N x e PA 



cos 17 vhx + sin 17 vhy - 

— sin 17 cos i vhx 

+ cos 17 cos i vhy 

+ sini v H z = e PE ■ v H 

, 

epE ■ 



epA ■ v H , 



(28) 



Vector e pa is directed from the Sun to the ascending node. 
epA ■ vh = Scosuj — /sinw is magnitude of vh compo- 
nent parallel with the line of nodes. The orbital plane is 
defined by its normal unit vector ejv. x epA = £pe- 
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Fig. 3. Secular time derivatives of S (see Eqs. 27) for a dust particle on prograde orbit in planar case. Origins of these 
Cartesian coordinate systems are in the Sun and vertical axes are aligned with the direction of the hydrogen gas velocity 
vector, ef/p and e^p are unit radial and transversal vectors in perihelion of the particle orbit (see text). 



£pE ■ vh — S'sinu; + Icosuj is magnitude of vh compo- 
nent perpendicular to the line of nodes and lying in the 
orbital plane. 



3. Theoretical discussion 

Eqs. (22)- (26) enable to deduce some properties of secular 
evolution of the dust particle under the action of the flow of 
interstellar gas. C = for a special case when the velocity 
of hydrogen gas vh lies in the orbital plane of the particle. 
In this case we get that the inclination and the longitude of 
ascending node are constant. Secular time derivatives for a 
special case of this kind (i = Cl = 0, vjjx = vhz = 0, vhy 
= vh in Eq. 27) have been derived in Klacka et al. (2009a). 
Secular time derivatives of a, e and ui are given in Klacka et 
al. (2009a) without generalization represented by Eqs. (27) 
and Fig. 2. 

Putting a = in Eqs. (22)-(26) one obtains a solution 
equivalent to the solution of Eq. (1) with the RHS inde- 
pendent of the particle's velocity (constant force). Eq. (22) 
yields constant semi-major axis for the unrealistic case a 
= 0. We will show, using Eq. (22), that secular semi-major 
axis decreases under the action of interstellar gas flow for a 
> 0. If secular increase of the semi-major axis would occur, 
then the value of curl braces in Eq. (22) would be negative. 
The value of 1 — s/1 — e 2 is always positive or zero. Thus, 
the curly braces could be negative only for I 2 — S 2 > 0. 
Since / and S are radial and transversal components of the 
constant vector vh in the perihelion (Fig. 2) of particle's 
orbit, we obtain maximal value of I 2 — S 2 for orbit orienta- 
tion characterized by S — 0. Using these assumptions, the 
minimal value (MV) in the curl braces is 



MV 




r2 i-Vi 



= l 



I 2 e 2 - 1 + VT^e 2 



H 
I 2 



T2 , 



l-VT 



> o . 



(29) 



The positiveness of MV means that the secular semi-major 
axis is a decreasing function of time. Eq. (22) was de- 
rived under the assumption that Solar gravity represents 
dominant acceleration in comparison to the interstellar gas 



disturbing acceleration. Since gravitational acceleration de- 
creases with square of a heliocentric distance, the assump- 
tion about the disturbing interstellar gas acceleration may 
not be valid far from the Sun. Behind the heliocentric dis- 
tance at which solar gravitational acceleration equals the in- 
terstellar gas acceleration, the particle's secular semi-major 
axis can also be an increasing function of time. The two ac- 
celerations are equal at heliocentric distance w 1 x 10 4 AU 
for dust particle with R = 1 /xm and g — 2 g/cm 3 , and, at 
heliocentric distance rj 7 x 10 3 AU for dust particle with 
R = 1 //m and g = 1 g/cm 3 . If Eq. (22) can be used, then 
secular time derivative of the semi-major axis is propor- 
tional to the value of semi-major axis (value of ^J-pj '/i a is 
independent of semi- major axis) . Our result differs from the 
Scherer's statement (semi-major of the particle can increase 
exponentially) at least due to an error in his expression for 
secular time derivative of magnitude of angular momentum. 

a 2 



Defining the function 10(e) = 1 — e" 
£ [0,1), present in Eq. (23), we obtain 



dw 
de 



/ 2 - vT 



: > 



VT 



(30) 



Thus, w(e) is an increasing function of eccentricity. We ob- 
tain that w(e) > for all e 6 [0,1), since w(e) is an in- 
creasing function of eccentricity and w(0) — 0. Hence, the 
sign of the second term in the square braces in Eq. (23) will 
depend on the sign of I 2 — S 2 . 

For small values of the eccentricity we get that the secu- 
lar time derivative of eccentricity is, approximately, propor- 
tional only to the first term multiplied by / in the square 
braces in Eq. (23). 

Eq. (24) yields that the argument of perihelion is con- 
stant for the planar case C = and for the orbit orientation 
characterized by S = 0. The longitude of the ascending node 
and the inclination are constant in this case. Hence differ- 
entiation of the second equation of Eqs. (27) with respect 
to time gives 



dl dui 
~dl~ ~ ~dt 



(31) 



Thus, / is also constant in this case. Since a is a small 
number, the first term in square braces in Eq. (23) is the 
dominant term on evolution eccentricity. 
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Now, let us consider the planar case (C = 0) with S ^ into averaged time derivatives of the quantities S, I and C 
0. The differentiation of the first of Eqs. (27) with respect defined by Eqs. (27), one obtains 
to time gives 



f=I^. (32) 
dt dt v ' 

This quantity can be averaged using Eq. (21). We get 



dS\ _j/du\ 
~dt/ ~ \~~dt 



(33) 



If g is a small number and I is not close to zero, then we 
can use the following approximations of Eqs. (23)-(24) 



de\ _ 3 cd -f H v H 
dt , 



doj 
~dt 



3 c D 7 g v H Fp S 
2 ffi e 



(34) 



(35) 



Inserting Eqs. (34)- (35) into Eq. (33) one can obtain 



dS_ 
~dt 



This equation leads to the differential equation 
de _ dS_ 

T ~~ ~s ' 

with the solution 
D 



(36) 



(37) 



(38) 



I? is a constant which can be determined from initial con- 
ditions. Thus, eccentricity is close to its minimal value if 
the major axis of the orbit is aligned with the direction of 
the hydrogen gas velocity vector. 

If we consider the system Eqs. (23)-(26) with a = 0, 

then we obtain (dS/dt) = — 3 Cd ^h vh vWm IS/ (2e) 
from the first of Eqs. (27). If we use Eq. (23) with a = 
0, then we get equation (de/dt) = — e (dS/dt) / S. This 
equation leads to the solution e — D / \S\ also for this non- 
planar case. But approximation a w is not allowed in Eqs. 
(24)- (26) since the first terms in square braces multiplied 
by C are multiplied by e, sinw and cos lj which can be close 
to zero. 

The case a — (the RHS of Eq. 1 is independent of 
the particle's velocity - constant force) is also significant 
for a motion of a rocket in a central gravitational field with 
a constant-reaction acceleration vector. We obtain the fol- 
lowing important result: if the constant acceleration can 
be considered as a perturbation acceleration to the central 
gravitation and S ^ 0, then equation e = D / \S\ holds 
during orbital motion of the rocket. This result fills up the 
results of Bctelsky (1964) and Kunitsyn (1966). 

Parameters S, I and C determine position of the orbit 
with respect to the hydrogen gas velocity vector. Therefore 
their time derivatives are useful for description of evolu- 
tion of the orbit position in space. Putting Eqs. (24)-(26) 



dS_ 

dt 



CD 1H VH S 





r 3/ 


a 









e 4 



e 4 -6e 2 +4-4(l-e 2 ) 3 / 2 



(dt) 




dC_ 
~dt 



CD lH V H 



+ ^ ( e 4 - 6e 2 + 4 



CD JH VH C Fp_ 
2 V jJL 



+ —(S 2 + I 2 ) 

V H 



3eC 2 
1-e- 

-4(1- 

3e7 

1-e 2 



3S 2 
H 

e 

o2^3/2 



al 

vh 



(39) 



(40) 



(41) 



Eq. (39) yields for the planar case (C = 0) 



dS_ 
~dt 



cp Jh vh S Fp_ 
2 V n 



3J 

e 



+ 



al 2 



2\3/2 



VHC 
CD 1H Vh S 



6e J + 4-4(l -e z ) 



*( - ^ + ^-fc(e) 



(42) 



Function fc(e) is a decreasing function of eccentricity for e 
G (0,1). This can be shown in a similar procedure as we 
used in Eq. (30) for the function w(e). The procedure must 
be used more than once for the function k{e). The function 
k(e) obtains values from lim^o k(e) = — 0.5 to lim e ^i k(e) 
= — 1, for e G (0,1). For the eccentricity e = 1 we obtain 
lim e _j.i a = oo. Hence we use maximal value of eccentricity 
e m for which Eq. (15) holds. The first term in parenthesis in 
Eq. (42) is minimal for e = e m for a given value of I. Value 
of the parenthesis in Eq. (42) will be maximally affected by 
the second term for e = e m . The second term in parenthesis 
yields for e — e m 



°I 2 w x vl 2 
k(e m ) > > 

v H v h 



\I\ >-\I\> 



(43) 



since a <C 1 and |/| < Vh- Thus, the sign of (dS/dt) de- 
pends only on the sign of the first term in the parenthesis 
in Eq. (42). Fig. 3 depicts unit vectors e_Rp = e_fj(/ = 0) 
and e-rp = er(/ = 0) in the orbital plane of a particle with 
prograde orbit. The hydrogen gas velocity vector vh for the 
planar case lies in the orbital plane. Direction of vh is also 
shown. The unit vector enp has direction and orientation 
from the Sun to perihelion of the particle orbit. In Fig. 3a 
epp lies in the first quadrant of the Cartesian coordinate 
system with origin in the Sun and vertical axis aligned with 
the direction of the hydrogen gas velocity vector. Both S 
and I are greater than 0, for these positions of the unit vec- 
tors epp and e^p- Hence (dS/dt) is negative. If epp lies 
in the second quadrant (Fig. 3b), then (dS/dt) is positive. 
If epp lies in the third quadrant (Fig. 3c), then (dS/dt) is 
negative. Finally, if epp lies in the fourth quadrant (Fig. 
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Fig. 4. Two evolutions of orbital elements of a dust particle with R = 2 pim and mass density g — 1 g/cm 3 under the 
action of interstellar gas flow. Evolution depicted by a solid line is calculated from the equation of motion. Evolution 
depicted by a dashed line corresponds to Eqs. (22)-(26). 



3d), then (dS/dt) is positive, as it is in the second quadrant. 
Therefore, the vector epp rotates counterclockwise in the 
first and the second quadrants, and, clockwise in the third 
and the fourth quadrants. If the vector epp is parallel with 
vertical axis in Fig. 3 (I = and (7 = 0), then Eq. (40) 
yields (dl/dt) > 0. Thus, positions of the vector epp par- 
allel with the vertical axis in Fig. 3 are not stable. However, 
if epp is parallel with the horizontal axis (S = and C = 
0), then (dS/dt) = and (dl/dt) = 0. Thus, e RP parallel 
with the horizontal axis yields stable positions of the vector 
epp . The stable position of the vector epp parallel with the 
horizontal axis and directed to the left in Fig. 3 is of the- 



oretical importance, only. In reality, no particles should be 
observed with perihelia in this direction. However, all unit 
vectors epp of particles in prograde orbits will approach 
the right direction in Fig. 3. For retrograde orbits we have 
to use transformation S — > S , I — > — / . We obtain (dS/dt) 
— > — (dS/dt). Hence, all unit vectors epp of particles in 
retrograde orbits will approach the left direction in Fig. 3. 
This result was obtained, in a different way, also by Scherer 
(2000). Scherer states that the approaches of unit vectors 
epp to one direction holds also for the non-planar case. 
However, the statement of Scherer is incorrect, in general. 
If / = and C ^ 0, then Eq. (39) implies that (dS/dt) is 
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proportional to SC 2 and not to SI. This leads to a behav- 
ior which differs from the behavior discussed above using 
Fig. 3. 

Scherer (2000, p. 332) furthermore states that orbital 
plane under the effect of the interstellar gas flow will be 
rotated into a plane coplanar to the flow vector vh, inde- 
pendent of initial position of the orbital plane. We showed, 
using numerical integrations, that the statement of Scherer 
is not correct. Retrograde orbits were not discussed by 
Scherer. We found an interesting orbit behavior. It depends 
on the orbit orientation with respect to the hydrogen gas 
velocity vector vh- If vh lies in a plane i = and vh is 
perpendicular to the vector e^p (in this case S = and 
(dS/dt) — 0), then the interstellar gas flow can change pro- 
grade orbit into a retrograde one (even more times for one 
particle) and inclination does not approach the value i = 
0. 



4. Numerical results 

4.1. Comparison of the numerical solution of equation of 
motion and the solution of Eqs. (22)-(26) 

We numerically solved Eq. (3) and the system of differen- 
tial equations represented by Eqs. (22)- (26). Solutions are 
compared in Fig. 4. We assumed that the direction of the 
interstellar gas velocity vector is identical to the direction 
of the velocity of the interstellar dust particles entering the 
Solar System. The interstellar dust particles enter the Solar 
System with a speed of about Voq = 26 km/s (Landgraf et 
al. 1999) and they are arriving from direction of A ec ; = 259° 
(heliocentric ecliptic longitude) and f3 ec i — 8° (heliocentric 
ecliptic latitude) (Landgraf 2000). Thus, components of the 
velocity in the ecliptic coordinates with x-axis aligned to- 
ward the equinox are vh — — 26 km/s (cos(259°) cos(8°), 
sin(259°) cos(8°), sin(8°)). As an initial conditions for a 
dust particle with R = 2 fxm and mass density g = 1 g/cm 3 
we used a in = 60 AU, e in = 0.5, uj in = 90°, Q in = 90°, i in = 
20° for Eqs. (22)-(26). The initial true anomaly of the dust 
particle was /j n = 180° for Eq. (3). Fig 4. shows that the 
obtained evolutions are in a good agreement. Evolutions be- 
gin separate as the eccentricity approaches 1. This is caused 
by the fact that approximation a -C 1, see Eq. (15), does 
not hold for large eccentricities. Detailed numerical solution 
of the equation of motion (Eq. 3) yields that the secular 
semi-major axis is a decreasing function of time also when 
eccentricity approaches 1. 

We can summarize, on the basis of the previous para- 
graph. If there is no other force, besides solar gravity and 
the flux of interstellar gas, then the semi-major axis of an 
interplanetary dust particle decreases and the particle can 
hit the Sun. However, the particle can hit the Sun also 
by another possibility: particle's eccentricity increases to 1. 
These mathematical possibilities probably do not occur in 
reality, since other forces can act on the dust particle and 
the interstellar gas is ionized below the heliocentric distance 
of about 4 AU. 

Let us return, once again, to the planar case (C = 0) in 
which 5 = and the dominant term in the square braces 
in Eq. (23), the term (3/2) /, is negative. Numerical inte- 
gration of Eq. (3) shows that if the eccentricity decreases 
to 0, then the argument of perihelion uj "shifts" its value to 
the value uj + n. This means that the negative value of / 




0.4 - 



2 4 6 8 10 
t[10 5 years] 




30 - 



2 4 6 8 10 
t [ 10 5 years ] 

Fig. 5. Orbital evolution of dust particle with R = 2 /xm, 
mass density g = 1 g/cm 3 and Q' pr = 0.75 under the action 
of the Poynting-Robertson effect (PR), the radial solar wind 
(SW) and the flow of interstellar gas (IG). 



changes to positive and the eccentricity begins to increase 
with the same slope. 

The approximative solution represented by Eq. (38) is 
in a good agreement with the detailed numerical solution 
of Eqs. (22)-(26) for the planar case with S ^ 0. This holds 
for the whole time interval, also for / close to zero. Eq. 
(38) holds, approximately, also for the non-planar evolution 
depicted in Fig. 4. In this case i is close to zero, vh ~ 
vhy and « 90° at the eccentricity minimum. Eq. (38) 
gives e « D / \vh coscj\. The evolutionary minimum of 
eccentricity occurs for the case when to is close to 180°. 
This is in accordance with Eq. (38). 
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4.2. Comparing of influences of interstellar gas flow, 
Poynting-Robertson effect and radial solar wind on 
dynamics of dust particles 

We have considered only the effect of interstellar gas flow, 
up to now. In reality, some other non-gravitational effects 
play non- negligible role. Thus, we want to compare the ef- 
fect of the interstellar gas flow with the other effects influ- 
encing dynamics of dust grains in the Solar System. For this 
purpose we included the Poynting-Robertson effect (P-R ef- 
fect) and the radial solar wind into the equation of motion. 
The P-R effect is electromagnetic radiation pressure force 
acting on a spherical particle (Klacka 2004; arXiv:astro- 
ph/0807.2915; arXiv:astro-ph/0904.0368). Equation of mo- 
tion for the P-R effect, the effect of the radial solar wind 
and the effect of gas flow has the form (e.g. Klacka et al. 
2009b) 



dv 
~dt 



en 



- c D -f H \v - v H \ (v - v H ) , 



V 

en, + - 

c 



(44) 



where decrease of particle's mass (corpuscular sputtering) 
and higher orders in v/c are neglected, c is the speed of light 
in vacuum. Parameter /? is defined as the ratio of the radial 
component of the electromagnectic radiation pressure force 
and the gravitational force between the Sun and the particle 
in rest with respect to the Sun: 



3L Q Q' 



WttcGMqRq 



5.763 x 10" 



Qm 



R [m] q [kg/m 3 



■ (45) 



L Q is the solar luminosity, L Q = 3.842 x I0 26 W (Bahcall 
2002), Q' pr is the dimensionless efficiency factor for radia- 
tion pressure integrated over the solar spectrum and cal- 
culated for the radial direction (Q pr = I for a perfectly 
absorbing sphere) and q is mass density of the particle. r\ 
is the ratio of solar wind energy to electromagnetic solar 
energy, both radiated per unit of time 



4-7rr u 



© 



JV 

E 



(46) 



where u is the speed of the solar wind, u = 450 km/s, 
m,i and rii, i = 1 to JV, are masses and concentrations of 
the solar wind particles at a distance r from the Sun, rj 
= 0.38 for the Sun (Klacka et al. 2009b). Four numerical 
integrations of Eq. (44) are shown in Fig. 5. We used dust 
particle with R = 2 fim, mass density q = 1 g/cm 3 and Q' pr 
= 0.75. We took into account only planar case when the 
interstellar gas velocity lies in the orbital plane of the dust 
particle (C = 0), for the sake of simplicity. Initial position 
is r in = (0, —90 AU, 0) and initial velocity vector is v in 
= (2 km/s, 0, 0). Orbital evolution is given by evolution 
of orbital elements calculated for the case when a central 
acceleration is defined by the first Kcplerian term in Eq. 
(44), namely — /u(l — f3)en/r 2 . The evolution depicted by 
the dash-dotted line is for the P-R effect alone (jh = 0, 77 
= in Eq. 44). The evolution depicted by the dotted line 
is for the P-R effect with the flow of interstellar gas (rj = 
in Eq. 44). The evolution depicted by the dashed line is 



for the P-R effect and the radial solar wind (7^ = in Eq. 
44). Finally, the evolution depicted by the solid line holds 
for the case when all three effects act together. 

The evolution of semi-major axis depicted in Fig. 5 
shows that the flow of interstellar gas is more important 
than the radial solar wind, as for the effects on the dynam- 
ics of the dust particle. 

The secular eccentricity is always a decreasing function 
of time for the P-R effect and the radial solar wind (e.g., 
Wyatt & Whipple 1950, Klacka et al. 2009b). The growth 
in eccentricity depicted in Fig. 5 is due to the interstellar 
gas. Fast decrease of the semi-major axis in Fig. 5 may 
also be, at least partially, caused by the fact that higher 
eccentricities decrease the value of (da/dt)pR and the P-R 
effect becomes stronger. The secular evolution of eccentric- 
ity can be also an increasing function of time if the flow of 
interstellar gas is taken into account. We have 



dt 




Mi - P) 

*(l1-Sg)(l-e|) 



x I 



1-4 



(47) 



if also Eq. (23) is used. The subscript (3 denotes that the 
central acceleration — /3)en/r 2 is used for calcula- 

tion of the osculating orbital elements. We remind that the 
transformation fi — > fi(l — f3) has to be done on the RHS 
sides of Eqs. (20)-(26). If we use definition of the oscu- 
lating orbital elements, then the physical central accelera- 
tion is given by the gravitational acceleration of the Sun, 
—\ienjr 2 . In this case, the secular evolution of eccentricity 
is given by Eq. (103) in Klacka (2004), assuming that ep is 
calculated from Eq. (47). 

If optical properties of the dust particle are constant, 
then the secular time derivative of the argument of perihe- 
lion equals to zero for the P-R effect and the radial solar 
wind (Klacka et al. 2007, 2009b). If the flow of interstellar 
gas is included into the equation of motion, then, even in 
the planar case, the secular time derivative of the argument 
of perihelion may not be equal to zero, in general (see Fig. 
5). 

The evolutions of eccentricity and argument of perihe- 
lion shown in Fig. 5 are significantly affected by the flow of 
interstellar gas. 

4.3. Dust ring in the Edgeworth-Kuiper belt zone 

Real situation in the Edgeworth-Kuiper belt zone may be 
much more complicated than the situation discussed in 
Sees. 4.1 and 4.2. In particular, gravitation of planets may 
have an important influence on dynamics of dust in the 
zone. For this reason we included gravitation of four ma- 
jor planets into the final equation of motion. Observations 
from Helios 2 during its first solar mission in 1976 (Bruno 
et. al. 2003) show that the angle between the radial direc- 
tion and the direction of the solar wind velocity does not 
significantly depend on heliocentric distance. If the value of 



P. Pastor et al.: Orbital evolution under action of fast interstellar gas flow 



9 



t = years 



100 




t = 25 x 10 years 



100 



100 



H 
r— 1 





. » 


-50 




-m 



t = 50 x 10 years 



100 



100 




100 



t = 75 x 1 years 



100 



100 



100 



-100 



100 -100 

Fig. 6. Time evolution of ring of dust particles with R = 2 fim, g=l g/cm 3 and Q' pr = 0.75 in the zone of the Edgeworth- 
Kuiper belt. The ring becomes eccentric in less than 10 6 years due to the interstellar neutral gas. Orbits of the particles 
are shown with black color and orbits of the planets are shown in gray. 



this angle is approximately constant, then the non-radial 
solar wind can also have an important influence on dynam- 
ics of dust in outer parts of the Solar System. We took into 
account the non-radial solar wind with constant value of 
the angle. Influence of precession of the rotational axis of 
the Sun on the non-radial solar wind was also considered. 
Hence, equation of motion of the dust particle has the form 



dv 
~dt 



+ 



r 

P 7} fj, 



v-e R \ v 
1 e R 

c J c 



„ \v — u\ (v — u) 
Q' pr cu r 2 1 11 1 

- c D j H \v - v H \ (v - V H ) 
V GMi (r r) 



E 



GM, 



(48) 



where u is solar wind velocity vector, Mj are masses of the 
planets and r*j are position vectors of the planets with re- 
spect to the Sun. The non-radial solar wind velocity vector 
was calculated from equation 



k x e R . 

u — u e b cos e + — sin e 

\k x e R \ 



(49) 



where e is the angle between the radial direction and the 
direction of the solar wind velocity and k is unit vec- 
tor with direction/orientation corresponding to the direc- 
tion/orientation of solar rotation angular velocity vector. 
Vector k for a given time can be calculated from 

k = (sinf2 s sini s , — cosf2 s sini s , cosi s ) , 

i s = 7° 15', fl s = 73°40' + 50.25" (t[years] - 1850) . (50) 
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Fig. 7. Evolution of orbital elements of dust particles in the Edgeworth-Kuiper belt zone during 96 numerical solutions 
depicted in Fig. 6. Initial values of orbital elements are: dj„ = 60 AU, e in = 0.1, uj in e {0°, 45°, 90°, 270°, 315°}, fij„ 
€ {0°, 90°, 180°, 270°}, i in e {5°, 10°, 15°} and f in = 0°. Evolution during the first 750 000 years is influenced mainly 
by interstellar gas, and, later on, mainly by gravitation of planets (see text). 



While Eqs. (49)-(50) are consistent with Klacka (1994) and 
Abalakin (1986), the value of e, e = 2.9°, used in our nu- 
merical calculations, is in accordance with Bruno et al. 
(2003) and Klacka et al. (2007). The observed neutral hy- 
drogen gas velocity vector in the ecliptic coordinates with 
x-axis aligned toward the equinox is vr = — 26 km/s 
(cos(259°) cos(8°), sin(259°) cos(8°), sin(8 )). As for the 
initial conditions of dust particles we did not use random 
positions and velocities. We assumed that putative dust 
ring in the Edgeworth-Kuiper belt has approximate circu- 
lar shape and contains lot of particles with approximately 



equal optical properties. We assumed that the ring contains 
such large amount of particles that one can choose, approx- 
imately, a given value of semi-major axis and a given radius 
of the particles. Accelerations caused by the P-R effect, the 
solar wind and the interstellar neutral hydrogen gas are in- 
versely proportional to particle's radius and mass density. 
Therefore, a large particle is less influenced by the non- 
gravitational effects than a small particle of the same mass 
density. The evolution of the large particle under action of 
the non-gravitational effects is slower than the evolution of 
the small dust grain. We used uniformly distributed initial 
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values of the argument of perihelion and longitude of the 
ascending node. Furthermore we assumed that particles in 
the ring orbit prograde in low inclination orbits. We used 
particles with R = 2 /im, g = 1 g/cm 3 and Q' pr — 0.75. 
Exact initial values of orbital elements were ai n = 60 AU, 
e in = 0.1, w in G {0°, 45°, 90°, 270°, 315°}, fl in e {0°, 
90°, 180°, 270°}, i in e {5°, 10°, 15°} and f in = 0°. Hence, 
we obtained 8 x 4 x 3 = 96 individual orbits. Results of 
numerical solutions of Eq. (48) are depicted in Figs. 6 and 
7. Fig. 6 depicts evolution of the dust ring viewed from per- 
spective. Orbits of the planets are also shown. Time span 
between various pictures in Fig. 6 is 250 000 years. As we 
can see, the ring becomes more and more eccentric because 
of a fast increase of eccentricity caused by the interstellar 
gas flow (see also eccentricity evolution in Fig. 7). Perihelia 
of orbits are shifted in accordance with the behavior dis- 
cussed in Fig. 3. This is caused by the facts that influence 
of interstellar gas is dominant and the solved problem is 
almost coplanar. The term multiplied by C 2 in Eq. (39) 
does not have large influence on the first term in the curly 
braces in Eq. (39), in almost coplanar case. Therefore, the 
lines connecting the Sun with the perihelia of particles or- 
bits are approaching the direction perpendicular to the in- 
terstellar gas velocity vector. Time evolution of the orbital 
elements of the dust particles in the ring, considered in Fig. 
6, is depicted in Fig. 7. Due to the approach of the peri- 
helia to one direction, the evolutions of the argument of 
perihelion uj, beginning with a given initial value u)i n , are 
divided into four branches. Each of them corresponds to 
one initial value of the ascending node. If the time is less 
than 750 000 years, then: i) the concentration of the parti- 
cles in the ring is smallest in the direction (from the Sun) 
into which perihelia of the orbits are approaching, and, ii) 
the concentration of the particles is greatest in exactly op- 
posite direction. If the time is greater than 750 000 years, 
then the orbits of particles in the dust ring are getting close 
to the orbits of the planets due to the increase of particles 
eccentricities. The situation after 750 000 years can be seen 
in Fig. 7. The dust particles with R = 2 ^urn, q = 1 g/cm 3 
and Q' pr = 0.75 are characterized by the value j3 ~ 0.216 
(see Eq. 45). For this value of fj, one obtains a rts 57.7 AU 
for the location of the exterior mean motion 3/1 resonance 
with Neptune; it follows from a = a P (1 - ^) 1/3 (3/l) 2/3 , 
where ap is semi-major axis of Neptune. We can see, from 
the evolution of semi-major axis in Fig. 7, that the secular 
semi-major axis is a decreasing function of time during the 
first 750 000 years. Thus the semi-major axis can evolve 
from an initial value of 60 AU to the location close to the 
mean- motion 3/1 resonance. Particles are influenced both 
by the vicinity of Neptune orbit and the exterior mean mo- 
tion 3/1 resonance with Neptune. Evolution during the first 
750 000 years is influenced mainly by the neutral interstel- 
lar hydrogen gas and, later on, mainly by the gravitation 
of planets. 

Inclusion of the P-R effect, the non-radial solar wind 
and the interstellar gas into the equation of motion of the 
dust particle without planets can stabilize the particle's or- 
bit. The stabilization is characterized by stable values of 
orbital elements. This stabilization is discussed in Klacka 
et al. (2009a) for sine = 0.05. The process of stabilization 
requires about 1 x 10 8 years for the dust particle with (3 
= 0.01. This time is not very sensitive to the efficiency fac- 
tor for radiation pressure Q' r . However, for lower values of 



Q' pr the stabilization occur with larger probability, because 
stabilization effect of the non-radial solar wind is stronger 
(see Eq. 48). If also planets are considered in the equation 
of motion, then the stabilizing value of dust particle eccen- 
tricity is usually sufficiently high to get the particle close 
to one of the planets during a long time span. As a con- 
sequence, gravitation of the planet can change orbit of the 
particle and cancel the stabilization process. 

5. Conclusion 

We investigated orbital evolution of a dust grain under the 
action of interstellar gas flow. We presented secular time 
derivatives of the grain's orbital elements for arbitrary ori- 
entation of the orbit with respect to the velocity vector of 
the interstellar gas, which is a generalization of several re- 
sults presented in Klacka et al. (2009a). The secular time 
derivatives were derived using assumptions that the accel- 
eration caused by the interstellar gas flow is small in com- 
parison with gravitation of a central object (the Sun), ec- 
centricity of the orbit is not close to 1 and the speed of 
the dust particle is small in comparison with the speed 
of the interstellar gas. These assumptions lead to secular 
decrease of semi-major axis a of the particle. The secular 
time derivative of the semi-major axis is negative and pro- 
portional to a. This result is not in accordance with Scherer 
(2000) who has stated that the semi-major of the particle 
increases exponentially. Scherer's statement is incorrect and 
our analytical result is confirmed by our detailed numerical 
integration of equation of motion (see also Fig. 4). 

If the hydrogen gas velocity vector vh lies in the parti- 
cle's orbital plane and the major axis of the orbit is not per- 
pendicular to vh, then the product of (secular) eccentricity 
and magnitude of the radial component of Vh measured in 
perihelion is, approximately, constant during orbital evolu- 
tion. 

We considered simultaneous action of the P-R effect, the 
radial solar wind and the interstellar gas flow. Numerical 
integrations showed that the action of the flow of inter- 
stellar gas can be more important than the action of the 
electromagnetic and corpuscular radiation of the Sun, as 
for the motion of dust particles orbiting the Sun in outer 
parts of the Solar System (see Fig. 5). Physical decrease of 
semi-major axis can be more than 2-times greater than the 
value produced by the P-R effect and radial solar wind. The 
evolution of eccentricity can also be an increasing function 
of time when we consider the P-R effect and the radial solar 
wind together with the flow of the neutral interstellar gas. 
This is also relevant difference from the action of the P-R 
effect and the radial solar wind when secular decrease of 
eccentricity occurs. Simultaneous action of all three effects 
yields that the secular time derivative of the argument of 
perihelion may not be equal to zero, in general. 

Gravitation of four major planets was also directly 
added into the equation of motion, see Eq. (48). This access 
correctly describes capture of dust grains into mean motion 
resonances with the planets. Our physical approach differs 
from the Scherer's approach (Scherer 2000), who has used 
some kind of secular access to gravitational influence of the 
planets. 

Assumption on an existence of dust ring in the zone of 
the Edgeworth-Kuiper belt is in contradiction with rapid 
increase of eccentricity of the ring due to an acceleration 
caused by the interstellar gas flow. Speed of the eccentricity 



12 



P. Pastor et al.: Orbital evolution under action of fast interstellar gas flow 



increase (time derivative of eccentricity) is roughly inversely 
proportional to the particle's size and mass density. As the 
eccentricity of the particles increases, the particles approach 
the planets. The particles in the ring are under the gravi- 
tational influence of the planets. The particles evolve also 
in semi-major axis and they can be temporarily captured 
into a mean motion resonance. The particles can remain in 
chaotic orbits between orbits of the planets, or, the particles 
are ejected to high eccentric orbits due to close encounters 
with one of the planets. Only particles with greater size 
and mass density should survive in the dust ring for a long 
time. 

A relevant result of the paper is that equation of motion 
in the form of Eq. (44) (or, Eq. 48) and Eqs. (45)- (46) have 
to be used in modeling of orbital evolution of dust grains 
in the Solar System. The influence of the fast interstellar 
neutral gas flow might not be ignored in general investi- 
gations on evolution of dust particles in the zone of the 
Edgeworth-Kuiper belt. 
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